Self-consistent model of the electron-positron pair production in the collision of laser 

pulses 



O 

(N 



Alexander Tarasenko 

Theoretical Physics Department, Belarusian State University, Belarus, Minsk, Nezavisimosti av., 4, 22003^\ 

(Dated: January 14, 2013) 

Pair production due to the Schwinger mechanism is modeled in a self-consistent approach taking 
into account the backreaction of the created particles on the field. The model uses the Boltzmann 
kinetic equation for the electrons and positrons and the Maxwell equations for the electromagnetic 
field, which are solved numerically in a certain approximation for the case of the linearly polarized 
plane wave collision. The results indicate three main stages of the collision process: the acceleration 
of particles, the deceleration of particles, and plasma oscillations. A linear equation describing the 
oscillations is obtained, the oscillation frequency is related to the frequency of the quasinormal mode 
of the plasma with the least damping. It is also shown that the oscillation frequency is of the order 
of the plasma frequency. 
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I. INTRODUCTION 

The electron-positron pair production in strong elec- 
tromagnetic field was predicted by Klein and Sauter (see 
Ref. [l| for a review), and modern laser facilities such as 
Extreme Light Infrastructure (ELI) 0| and High Power 
laser Energy Research facility (HiPER) |3| offer a pos- 
sibility to check these predictions. Potential intensities 
that can be reached at these facilities are up to 10 26 
W/cm 2 IS], which allows to observe physical processes 
taking place in an electric field close to the critical quan- 
tum field E c = m^c 3 /eH, which corresponds to the inten- 
sity of the order of 10 29 W/cm 2 . 

The phenomena taking place in a superstrong laser 
field have attracted much attention of researchers re- 
cently. Various processes have been investigated such 
as the electron-positron-photon cascade formation (U - !!] 
and different quantum electrodynamical reactions includ- 
ing the nonlinear Breit- Wheeler and Bethe-Heitler effects 
, photon- photon scattering Q , and muon pair produc- 
tion Q. Fedotov et al. Q emphasize that if a cascade 
develops the influence of the particles on the initial elec- 
tromagnetic field has to be taken into account. Such self- 
consistent model of cascade development in the field of 
colliding linearly polarized laser pulses has been recently 
proposed by Nerush et al. @. 

A self-consistent approach to pair production has been 
developed in a simpler model treating the electric field 
as uniform and considering only the Schwinger mecha- 
nism of pair production. Both the model based on the 
semiclassical quantum field theory and the one based 
on the the kinetic description of the electron-positron 
plasma (Tll - [T3| predict plasma oscillations [l[ . In a suffi- 
ciently large time when the pair production almost stops 
these oscillations can be interpreted as natural oscilla- 
tions of the plasma. We may expect that the plasma 
formed during the collision of laser pulses also has nat- 
ural oscillations that show themselves after the collision. 



The goal of the present paper is to model the evolution 
of the electron-positron plasma and the electromagnetic 
field in the collision of two linearly polarized laser pulses, 
taking into account the backreaction of the particles on 
the field. 

The system (the electron-positron plasma and the 
field) will be modeled under the following assumptions: 
each of the pulses is a plane linearly polarized electro- 
magnetic wave. The pairs are produced only due to the 
Schwinger mechanism, i. e. we neglect the radiation reac- 
tion force acting on each particle and thus the possibility 
of cascade formation. The electric field is of the order of 
or less than the critical one, and the characteristic scale 
of its variation is much less than the Compton wave- 
length of the electron A e = h/m e c, which implies that 
the particles are produced almost at rest and the pair 
production rate is described by the Schwunger formula. 
The electron-positron plasma is considered to be ideal, 
nondegenerate and collisionless. The plasma is treated 
in the ultrarelativistic approximation, i. e. we assume 
that in a negligible time after the creation the particle 
moves at the speed of light. 



II. KINETIC DESCRIPTION OF THE 
ELECTRON-POSITRON PLASMA 

Consider an electromagnetic field with the potentials 
(p = 0, A = A y (t,x)e y , where (x,y,z) are the Cartesian 
coordinates. The electric field strength E and the mag- 
netic induction B are directed along y-axis and z-axis, 
respectively. The only nonzero components of these vec- 
tors are 



F IdAy 



dAy 

dx 



(1) 



Consider the fields for which A y (t,—x) = A y (t,x) and 



thus E y (t,—x) = E y (t,x) and B z (t,- 



-B z (t,x). 
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Such field can model the collision of waves at the point 
x = where the magnetic field vanishes. We will consider 



2 



the following initial field 



erf 



7T / . / ct — x — a 

— E a\ erf = 

2 I V V2a 



ct + x — a 



(2) 



V2o 



where erf x is the error function 



erf (x) = —= 

\/7T 



dt. 



(3) 



Field |2| is a linear combination of two gaussian 
wavepackets with the amplitude (the maximum value of 
the electric field strength) Eq and the characteristic size 
a. The wavepackets are initially located at the distance 
2a and collide at the moment t = a/c, the characteris- 
tic duration of the collision being a/c. The maximum 
intensity of each wave is 



cE 2 

W = — - = 4.65 x 10 29 

47T 



W/cm 2 



(4) 



In order to ensure the applicability of the Schwinger for- 
mula for pair production and the ultrarelativistic approx- 
imation we assume that the following inequalities are sat- 
isfied: 



E <E C , <7>A e , e£ , CT>m e c 2 



(5) 



Denote f e (t, x, p) and f p (t, x, p) the distribution func- 
tions of electrons and positrons, respectively. Since the 
particles are created at rest they move in the xy plane in 
field ([TJ. The Boltzmann kinetic equation for electrons 
reads 



dfe , df e 



dt 



dx 



-v y B z 



dp x 



< \E V ~ -^B^J ^ = S{E y ,B z )6{p), 



(6) 



where v is the electron velocity, S(E, B) is the source de- 
scribing the pair production due to the Schwinger mech- 
anism fl4f: 



S(E,B) 



S(E) 



dN 
dtdV 
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9(E 2 - B 2 )S(VE 2 - B 2 ), 



47T 3 fr 4 \E, 



E 



exp 



nEc 
E 



(7) 



and 6(x) is the Heaviside theta- function. Note that in 
general the pair production rate depends on two field 
invariants: E 2 — B 2 and EB. For field ((TJ) the sec- 
ond invariant is identically zero, thus the pair produc- 
tion is described by the same expression as the one in 
the electric field E' = V ' E 2 — B 2 . The theta-function in 
expression ([7]) prohibits pair production in a field that 



becomes a purely magnetic one in a certain reference 
frame. The source in equation ^ is proportional to 
the delta-function of the particle momentum which is 
a consequence of the assumption that the particles are 
produced at rest. It is valid if the characteristic energy 
obtained by a particle during the acceleration in the field 
is much greater than its initial energy, which is guaran- 
teed by inequalities ([S]). 

The distribution function depends on five variables, 
which complicates the numerical solution of the kinetic 
equation. Therefore we replace the exact equation by an 
approximate one neglecting the terms proportional to v x 
and B z in equations ([6]) . This approximation can be ex- 
plained in the following way: the pair production rate 
in a subcritical field has a sharp maximum around the 
maximum point of the electric field, i. e. x = 0. The 
magnetic field in this region is small due to the symmetry 
relation (B z (t, — x) = —B z (t,x)), and the velocity com- 
ponent v x is appears only due to the magnetic field and 
thus is also small as compared to the component v y . We 
note that these conditions do not guarantee that the par- 
ticle remains in the region with weak magnetic field all 
the time, and the particle really leaves this region. How- 
ever, we will consider only the approximate model here 
assuming that it correctly predicts the qualitative prop- 
erties of the system (see also Discussion section). The 
substitution f e (t,x,p) = F e (t,x,p y )5(p x )6(p z ) reduces 
the approximate equation to 



dF e 
dt 



eE. 



dFe 
dp v 



S(E y ,B z )S( Pv ). 



(8) 



There is no need in solving the kinetic equation for the 
positrons since f p (t,x,p x ,p y ,p z ) = f e (t,x,p x ,-p y ,p z ). 
The Maxwell equations reduce to 



<9 2 A, 1 d 2 A„ 



dx 2 



dt 2 



4-7T . 



(9) 



where the only nonzero component of the conduction cur- 
rent is 



j y (t,x) 



-Ice 



sgap y F e (t,x,py)dpy, (10) 



where sgnx is the signum function. There is no polar- 
ization current on the right side of equation (|9]) since it 
is identically zero for source ^ in the ultrarelativistic 
approximation. 

We solved equation system (|5))- ([TU)) numerically with 
the following initial and boundary conditions 

A y (0,x) = AfM, §- t A y (0,x) = ^4 0) (0,x), 

A y (t.±X ) = 0, F e {0,x,p v )=0, 
F e (t,±X ,p y ) = 0, F e (t,x, ±p max ) = 0, 

(11) 

where Xq and p m ax are sufficiently large numbers deter- 
mining the boundary of the integration domain. We fixed 
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the parameters a — 5(7 and a = 10~ 4 cm (which corre- 
sponds to the laser wavelength of the order of 1/im), and 
varied Eq from 0.065i? c to 0.23E C . A system of units 
with a = c = E c = 1 was used in calculations, so that 

A e = 3.86 x 1CT 7 , e = aE c X 2 e = 1.09 x 10" 15 , (12) 

where a — e 2 / he is the fine structure constant. The 
corresponding intensities of the initial waves Q var- 
ied from 1.96 x 10 27 W/cm 2 to 2.46 x 10 28 W/cm 2 . The 
first value is approximately the pair production thresh- 
old: for Eo = E m in = 0.0644_E C the pair production 
rate is S(E m i n ) = ca~ A . The second value, as the nu- 
merical solution showed, corresponds to a regime when 
the backreaction of the particles significantly changes 
the initial field. Note that the adiabaticity parameter 
r\ = aE n /X e E c [l2| is large in the range of parameters 
i?o and a used in calculations. The condition r\ ^> 1 
is necessary to use the Schwinger formula for the pair 
production rate. 

Equation system (|5|)- (ITTJ1) was solved in the region < 



t < T , -X < x < X 



< Py < Pn 



where 



To = -^o = 15, Pmax = 0.5. The finite difference method 
was used, the grid having N x = 251 points in the x- 
direction and N p = 401 point in the p^-direction. Due 
to strong dependance of the pair production rate on the 
field strength it was sufficient to solve equation ([8j in 
the range —HAx ^ x ^ llAir, where Ax = 0.12 is the 
grid spacing. Since the source in equation ([8]) is singular 
(proportional to the delta-function) the equation was not 
solved directly. The solution was expressed as 



F e (t,x,p y ) 



)(- Py )F^(t,x,p y ) 

'(p y )F^(t,X, Py ), 



(13) 



and both functions Fe (t,x,p v ) were obtained by solv- 
ing the corresponding homogeneous equation. Equation 
(|8j itself was used for solution joining: 



F^-\t,x,0)-F^(t,x,0) = -±-S{E y ,B 



eE v 



(14) 



The electric field E y at different moments of time for 
a specific solution with Eq = 0.15 is shown in Fig. [T] 
It follows from the figure that the backreaction of the 
particles on the field significantly changes the latter. It 
is also evident that at some points the electric field has 
the opposite direction to the direction of the initial field 

©• 

Equation system (|5j)- (fTU|) is self-consistent, which leads 
to the conservation of the total energy £ = £f+£ p , where 



1 

8~~- 



(E 2 y +B 2 z )dx 



(15) 



is the field energy and 



£p = 2 J J \p y \F e dp y dx 



(16) 
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1. (color online) Electric field E y (t,x) as a function of 
t = (a), t = 0.5T (b), and t = T (c) for the solution 
E = 0.15 (black) and for initial field © (red). 



is the energy of particles, both related to a unit area in 
the yz plane. Fig. [5h. shows the energy of particles as a 
function of time. It follows from the figure that at first a 
significant part of the total energy (about 90% for E > 
0.1) is transmitted to the particles. Then the electric 
field changes its direction and a large part of the energy 
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of particles is returned back to the field. For Eq > 0.1 
the particles possess from 10% to 40% of the total energy 
in the final state, the percentage increasing with Eq. The 
calculation shows that the total energy is conserved with 
the maximum deviation not exceeding 2% for all values 
of E . 

The concentration of electrons (or positrons) is given 

by 



= n p = / F e dp y . 



(17) 



The numerically calculated electron concentration is pre- 
sented in Fig. [2b- The figure shows that the pairs are 
produced in a region with a characteristic size a centered 
at the point x = 0. Note that plasma with a character- 
istic Lorenz factor of particles equal to 10 5 and a con- 
centration of the order of 10 25 cm -3 can be considered as 
ideal and nondegenerate, which was initially assumed. 
The total number of pairs per unit area can be written 

as 



N 



n e dx 



F e dp y dx. 



(18) 



The evolution of the total number of pairs is shown in 
Fig. |3] It follows from the figure that pairs are produced 
in significant quantities only at certain moments of time 
which, as the analysis shows, approximately coincide 
with the maximum points of the electric field strength 
\E y \ at the origin x = 0. 

We also found the number and energy of pairs in a 
non-self-consistent model the basic equations of which 
consist of equation (jHJ and ((9]) without the source on 
the right side. The self-consistent model always predicts 
smaller values both for the number of particles and their 
energy, the difference in predictions increasing with Eq: 
for Eq = 0.065 the ratio of the numbers and energies of 
particles in the final state in self-consistent and non-self- 
consistent models are 0.99 and 0.94, respectively, and for 
E a = 0.23 these ratios are 6.0 x 10~ 7 and 8.3 x 10~ 8 . 



III. LINEARIZED MODEL 

Numerical solution gives the evidence of plasma oscil- 
lations after the collision, which are present in Fig. []}:. 
When the magnitude of the electric field is sufficiently 
small we expect the oscillations to be linear. In this case 
the quantity A y (t, x±) as a function of time at a fixed 
point x — X\ has to be of the following form: 

A y (t, Xl ) = Ay + Aqc-^-^ sin (« (t - t a )), (19) 

where the parameters A y , Aq, ujq, t a , A do not depend on 
time. The potential (Ti"9|) possesses two properties: firstly 
the moments of time t n form an arithmetical progres- 
sion, where t n are the roots of equation E y (t n ,x\) = 



enumerated in ascending order; secondly the quantities 
(A y (t n+ i, x\) — A y (t n , xi)) form a geometric progression. 
We chose xx = and found that these two properties 
hold for all values of Eq in the range 0.065 < Eq < 0.23 
if we omit several first moments t n at which the oscilla- 
tions are still nonlinear. The common difference and the 
common ratio of the progressions allow to calculate the 
oscillation frequency loq and the parameter A. Particu- 
larly, for E = 0.15 



ujq = 3.49, A = 0.785. 



(20) 



The linear regime of oscillations is described by linear 
equations that follow from equation system (|8|)-(rT7J|). In 
this regime the pair production is negligible, therefore 
equation © has a solution 

F e (t,X,Py) = * (x, Py - -Ay(t,X)J , (21) 



where ^(x^py) is an arbitrary function. Equation (|21j) 
implies that the plot of F e versus p y for a fixed value 
of x shifts along p y -axis with time, its form remaining 
constant. Introduce the following notations: 



Ay(t, X) = Ay(x) + 5Ay(t, x) , 
Ay[x) — Ay(00, X). 

In the linear in SA y approximation 



(22) 



j y (t,x) = ]y(x) -4e 2 5A y {t,x)t> (x,-^A y (x)j , (23) 



where j y (x) — j y (oo,x). Substitution of equation ([23 
into © yields 



d 2 1 d 2 



dx 2 c 2 dt 2 



V(x) SA y = 0, (24) 



V(x) 



Wire* 



F e (oo,x,0), 



(25) 



where relation (|2"Tj) and the fact that SA y = is a solution 
of the linearized equation were used. 

Equation (|24p is a wave equation with the potential 
V{x). Numerical solution of equation system (jg])- (Ti"U)) 
together with equation (|2"5|) allows to find the poten- 
tial for arbitrary values of the parameters. The analysis 
shows that the potential is always positive, symmetric 
(V(—x) = V(x)) and tends to zero at the infinity. For 
Eq < Ei ~ 0.15 the potential has a single maximum at 
the point x — 0. If Eq > Ei then the maxima of the 
potential are located symmetrically with respect to the 
origin and x = is the minimum point. Equation (|24j) de- 
scribes the penetration of a wave through a potential bar- 
rier. The asymptotic form of the function SA y for large 
t can be described in terms of the quasinormal modes 
of oscillations [l5[ , which are analogous to the quasidis- 
crete states of a quantum system. For sufficiently large 
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FIG. 2. (a) Energy of particles as a function of time for the solution with Eo = 0.15 obtained from equation (I16[) (solid line) 
and as £/(0) — £f(t) using equation (|15[) (shown by points). The coincidence of both plots implies the conservation of the total 
energy. The latter is equal to £ = 1.24 x 10 21 erg x cm" 2 , (b) The concentration of electrons as a function of x for the same 
solution. 
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FIG. 3. Number of pairs N as a function of time for the 
solution with Eo = 0.15. 



t the quantity 8A y is approximately proportional to the 
quasinormal mode with the least damping. The quasi- 
normal mode is a solution of equation (fM]) of the form 
5 A y = ip(x) exp (—iuit), where uj = ujq — iX (the param- 
eters ujq and A have the same meaning as in expression 
(fTTHO . The function if) satisfies the following equation 



V(x) ib(x) = 0, 



(26) 



where a prime denotes differentiation of a function with 
respect to its argument, together with the boundary con- 
ditions 



exp 



H*>). 



ip'i-Xx) =-%- exp ( -i-Xi 
c V c 



(27) 



where X\ is sufficiently large so that the inequality 
V(±Xi) <C uj 2 /c 2 holds (for an exact quasinormal mode 
X\ = oo). We found the quasinormal mode frequency 
numerically integrating equation (|26[) with initial condi- 
tions (|27l) for a fixed value of uj, and then solving equa- 
tion (|28|) as an algebraic one. The initial condition re- 
quired for the numerical solution of equation (|28[) was 
w ^ = cy/V rnax: where V max is the maximum value of 
the potential (the first quasidiscrete state is located near 
the maximum point of the potential energy, see also [15j ) . 

The calculated parameters loq and A are presented in 
Fig. 01 The data shown in the figure were obtained 
similarly to the values (|2U)) . i. e. by solving nonlinear 
equation system ([8|)- fll"0l) . The values of w and A were 
calculated only for Eq > 0.1 since in the opposite case 
we could not observe enough oscillations. The results 
of the quasinormal mode calculation lead to the same 
values of the parameters with a maximum deviation of 
less than 2% both for ojq and A. The plots (Fig. 2]) show 
a nontrivial dependance of the parameters on the field 
strength, particularly, no asymptotic behavior both for 
large and for small values of E can be revealed. The 
same conclusion is expected about the dependance of ujq 
and A on a. 

The oscillation frequency ujq can be estimated for 
Eq < E\ given the characteristic Lorenz factor 7 and 
the concentration of electrons n e . Using the characteris- 
tic value of the momentum m e cy and equation (|17p we 
obtain: n e (t,x) <~ F e (t,x,0)m e cry. The approximation 
ojq ~ Cy/V(Q) (for Eq < Ei the maximum value of the 
potential is V(0)) and formula (|25]l lead to the estimate 



w 



167re 2 n e 
m e 7 



(29) 



- i-il>{Xi) = 0, 
c 



Benedetti et al. [13( introduce the notion of the relativis- 

(28) 

tic plasma frequency uj p , which is y/2 times less than the 
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FIG. 4. Oscillation frequency ojq (a.) and parameter A (b) as a function of Eq. 



right side of equation ([25)) (it should be noted that in 
[13| all electrons are assumed to have the same Lorenz 
factor, and 7 is understood as its value). Estimate ([25)1 
shows that the oscillation frequency is of the order of 
the plasma frequency. The parameter A is related to the 
quantity V"(0) [H[ and thus is more difficult to estimate. 



IV. DISCUSSION 

The results of the modeling suggest the following sce- 
nario of the laser pulse collision. Pair production begins 
approximately at the moment of collision and occurs in 
a small region the size of which is of the order of the 
laser wavelength. The created particles are accelerated 
to ultrarelativistic energies by the electric field. The par- 
ticles produce electric current which causes the electric 
field to change its direction, and the energy of particles 
is returned back to the field. This process can repeat 
several times. After the magnitude of the field becomes 
sufficiently small the particle production almost stops, 
the field energy remains approximately constant, and the 
plasma reaches the equilibrium state through oscillations. 

The advantage of the considered model is that the ki- 
netic equation for the electrons is reduced to an equation 
in three variables, and even this equation has to be solved 
in a small region close to the origin (|x| < a). However, 
this allows to obtain only qualitative but not quantita- 
tive description of the system. The direct solution of 
the equation of motion shows that the electron leaves 
the region where the magnetic field is weak in time (the 
escape time is of the order of cr/c). This leads to the 
expansion of the plasma cloud, which was also obtained 
by Nerush et al. [6]. Moreover, the radiation losses of a 
particle become important since the particles are accel- 
erated to ultrarelativistic energies in subcritical fields (in 
agreement with [l2j )■ A self-consistent description of ra- 
diation requires taking into account the existence of pho- 



tons and the possibility of pair production by them. The 
processes of radiation and pair production by a photon 
result in the possibility of the electron-positron-photon 
cascade formation. Nevertheless, we expect that these 
two effects do not change the general conclusions drawn 
above, particularly about the existence of plasma oscilla- 
tions. The expansion of the plasma is expected to reduce 
the oscillations: the oscillation frequency has to decrease 
due to the decrease in particle concentration (f!?9")) , and 
the parameter A has to increase since the plasma cloud 
loses energy not only due to the radiation but also due 
to the particle loss. The cascade formation is expected 
to increase the oscillation frequency as it increases the 
particle concentration (it is more difficult to draw a con- 
clusion about the parameter A: it has to decrease if the 
maximum of the potential V(0) becomes less sharp and 
thus V"(0) decreases, but it is not clear whether such 
effect takes place). 

A model that is free of these two approximations but 
does not take into account the Schwinger mechanism of 
pair production was developed by Nerush et al. [(J . Pa- 
per [6| considers the cascade formation in the collision 
of linearly polarized laser pulses having finite width in x 
and y-directions. The model is based on the combination 
of the particle-in-cell (PIC) and Monte-Carlo algorithms, 
and is self-consistent. However, neither energy transmis- 
sion from the particles to the field nor plasma oscillations 
were obtained as the result of modeling. We can explain 
it by the fact that the numerical solution was stopped 
when the energy of all particles was roughly 40% of the 
total energy. The results of the present paper suggest 
that the energy transmission to the field starts when the 
energy of particles becomes approximately 90% of the 
total energy, and the oscillations start even later. 

Plasma oscillations were obtained in the uniform 
plasma model [HI, [l3j]. The plasma oscillations in this 
model in the weak field limit are interpreted as the lon- 
gitudinal Langmuir waves having an infinite wavelength 
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and thus the frequency equal to the plasma one. The elec- 
tromagnetic radiation is impossible in this model thus the 
field energy periodically vanishes being totally transmit- 
ted to the particles, and the oscillation damping becomes 
negligible for large t. In the given paper the geometry 
of the field is different thus several conclusions change: 
the field energy is not totally transmitted to the parti- 
cles, the percentage of the transmitted energy depending 
on the intensity of the initial waves; the oscillations are 
damped due to the electromagnetic radiation. However, 
both models predict the same estimates of the param- 
eters: the oscillation frequency is of the order of the 
plasma frequency (l3j |. the characteristic Lorenz factor 
of the particles is of the order of 10 5 for Eo ~ 0.15 fl2j |. 



V. CONCLUSIONS 

In the given work a self-consistent model of pair pro- 
duction in the collision of two linearly polarized plane 
electromagnetic waves is presented. The general scenario 
of collision is described; it consists of three stages: the ac- 
celeration of particles, the deceleration of particles, and 
plasma oscillations. The linearized model describing the 
plasma oscillations is obtained from the general nonlin- 
ear one. The parameters of oscillations are related to 
the complex frequency of the quasinormal mode with the 
least damping. Estimate of the oscillation frequency ([29]) 
is obtained. The effects that are not taken into account 
in the model are considered, and their role is discussed. 
The results of this paper are compared with the results 
obtained using different models. 
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